Hydrodynamic Correlation Functions of a Driven Granular Fluid in Steady State 
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We study a homogeneously driven granular fluid of hard spheres at intermediate volume frac- 
tions and focus on time-delayed correlation functions in the stationary state. Inelastic collisions are 
modeled by incomplete normal restitution, allowing for efficient simulations with an event-driven 
algorithm. The incoherent scattering function, -Fi nco h(g,i), is seen to follow time-density superposi- 
tion with a relaxation time that increases significantly as volume fraction increases. The statistics 
of particle displacements is approximately Gaussian. For the coherent scattering function S(q, ui) 
we compare our results to the predictions of generalized fluctuating hydrodynamics which takes 
into account that temperature fluctuations decay either diffusively or with a finite relaxation rate, 
depending on wave number and inelasticity. For sufficiently small wave number q we observe sound 
waves in the coherent scattering function S(q, w) and the longitudinal current correlation function 
Ci(q,u)). We determine the speed of sound and the transport coefficients and compare them to the 
results of kinetic theory. 

PACS numbers: 61.20.Lc, 51.20.+d, 45.70.-n, 47.57.Gc 
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I INTRODUCTION 

The long wavelength, low frequency dynamics of gran- 
ular fluids is frequently described by phenomenological 
hydrodynamic equations [l-|4| . In contrast to a fluid com- 
posed of elastically colliding particles, the total energy 
of the system is not conserved, implying a finite decay 
rate of the temperature in the limit of long wavelength. 
Hence, strictly speaking, the temperature is not a hydro- 
dynamic variable. More generally, the scale separation 
required by hydrodynamics has been questioned [5j]. If 
the system is not driven, the homogeneous state is un- 
stable [6] and large spatial gradients develop — invalidat- 
ing a hydrodynamic approach. A third point of criticism 
refers to the pressure in the Navier-Stokes equation. Clo- 
sure of the hydrodynamic equations requires an equation 
of state to express the pressure [7f in terms of density 
and temperature. However, an equation of state is ex- 
pected to exist only in an equilibrium state. Given these 
problems, the hydrodynamic approach has been mainly 
restricted to small inelasticity, such that the decay rate of 
the temperature is small, the time for the build-up of spa- 
tial inhomogeneities is long, and an equation of state is 
approximately valid. In this limit, kinetic theory has pro- 
vided a basis for the hydrodynamic equations and give n 
explicit expressions for the transport coefficients |8h1(|. 

Driving a granular fluid allows for a compensation of 
the energy which is dissipated in collisions, such that a 
non-equilibrium stationary state (NESS) is reached. In 
experiment, the driving is frequently performed by shear- 

or by driving 



ing ll|, Il2j . with vibrating walls [1 



scales, which might be generated by driving through the 
boundaries, when the agitation decays over a character- 
istic length, e.g. the width of a shear band. Hence in the 
follo wing, w e consider a homogeneously driven granular 
fluid (18I422J and set out to investigate the validity of the 
hydrodynamic approach in the NESS. 

We present in this paper results for a homogeneously 
driven system of hard spheres with moderate inelasticity, 
parametrized by a coefficient of restitution e = 0.8, 0.9 
and 1.0 (clastic). We use event-driven simulations to 
focus on the dynamics of the system. Previous stud- 
ies of correlation functions for granular fluids have been 
either on density-density correlations at the same time, 
such as the structure factor S(q) [Hi |l9( and the pair 
correlation function [23[, or on velocity- velocity corre- 
lations [21| (and references therein). We focus here on 
time dependent spatial correlations [22| at volume frac- 
tion 0.05 < rj < 0.4 and compute the incoherent and 
coherent intermediate scattering functions. 

The former entails information about the motion of a 
tagged particle, which is expected to be diffusive at long 
times. We find that the incoherent scattering function is 
well approximated by a Gaussian and obeys time-density 
superposition. The divergence of the relaxation time as 
a function of r\ occurs not only for the elastic case but 
also for the inelastic case, consistent with the results of 
Reyes et al. [H| and J22|. 



the system homogeneously [lfj, |l_7| . To test a hydrody- 
namic approach in the NESS we want to avoid new length 



The coherent correlations reveal the collective dynam- 
ics of the fluid: damped sound waves and relaxation of 
temperature fluctuations. We determine the dynamic 
structure factor S(q,w) and compare our data quanti- 
tatively to the predictions of van Noije et al. [19( using 
fluctuating hydrodynamics. The agreement between sim- 



ulations and theory is quite good: damped sound waves 
are indeed observed for small wave numbers q and the 
velocity and damping of sound can be determined. Tem- 
perature fluctuations are found to decay either diffusively 
or with a finite rate, depending on q. The transport coef- 
ficients are compared to the predictions of kinetic theory 
and found to agree well. 

In the following we specify model and simulation de- 
tails in Sec. II. Subsequently, in Sec. Ill we discuss the 
incoherent scattering function, the mean square displace- 
ment, and the diffusion constant. Data for the intermedi- 
ate coherent scattering function and longitudinal current 
correlation function are presented in Sec. IV. A hydrody- 
namic model, which was first introduced in Ref. [18| . is 
discussed in Sec. V and compared to the simulation data 
for the coherent scattering function in Sec. VI. 



II MODEL AND SIMULATION DETAILS 

We investigate a system of N monodisperse hard 
spheres of diameter a and mass m at volume fraction 
r) = N Jy . The time evolution is governed by instanta- 
neous inelastic two-particle collisions. We consider here 
only the simplest model of an inelastic two-body colli- 
sion, described by incomplete normal restitution. The 
change of the relative velocity g := vi — V2 of the two 
colliding particles is given by 



(g • n)' = -e (g • n) , 



(1) 



where primed quantities indicate post-collisional veloci- 
ties and unprimed ones refer to precollisional ones. The 
unit vector n := (i"i— r^)/ |(i"i — r2)| connects the centers 
of the two spheres, and e =const. G [0, 1] denotes the co- 
efficient of normal restitution, with e — 1.0 in the elastic 
limit. The postcollisional velocities of the two colliding 
spheres are given by 



Vl 



v 2 



(1+f) 
2 

(1 + e) 



(n-g)n 
(n-g)n. 



(2) 
(3) 



Due to the inelastic nature of the collisions, we have 
to feed energy into the system in order to maintain a 
stationary state. The simplest bulk driving 24| consists 
of a kick of a given particle, say particle i, instantaneously 
at time t, which corresponds to 



kicking the particles randomly with amplitude «Dr and 
frequency / Dr . 

If a single particle is kicked at a particular instant, mo- 
mentum is not conserved. Due to the random direction 
of the kicks the time average will restore the conservation 
of global momentum, but only on average. Momentum 
conservation is known to be essential for the dynamic cor- 
relation functions in the limit of long wavelength and long 
times. Hence we choose a driving mechanism in which 
pairs of particles are kicked in opposite directions [25]. 
The pairs are fixed globally so that the total momentum 
is conserved at each instant of time. Denoting the part- 
ner of particle i by p(i), the random force correlation is 
given by 

& (a) (t)£ i °V)) = £oto,< - <W V(* - 0- (6) 

It is also possible to ensure momentum conservation on 
small scales by choosing pairs of neighboring particles 
and by kicking them in opposite directions. However this 
is not pursued here. 

For the event-driven simulations we use the optimized 
algorithm of Lubachevsky 26] adapted to granular me- 
dia [21]. To avoid the inelastic collapse we use the tech- 
nique of virtual hulls around the spheres as described 
in [2l|. Particles are colliding elastically when they 
are a diameter a apart and the dissipation takes place 
when the colliding spheres are receding and separated by 
(l + 10~ 4 ) a. With the appropriate choice of v^ T fui = £o 
we ensured constant temperature and in all following we 
chose units such that m = a = T = 1. All simulations 
were with a cubic box and periodic boundary conditions. 
The simulation results of Sec. Ill are for TV = 200000 and 
two independent simulation runs [31] , whereas for the re- 
sults of Sec. IV and Sec. VI we needed more statistics and 
therefore used N = 10000 and 100 independent simula- 
tion runs. In each set of simulations we first equilibrated 
at e = 1.0 at the desired volume fraction, followed by a 
relaxation to a stationary state at e ^ 1 (achieved with 
a simulation run of at least 100 time units) and consec- 
utive production runs. Independent configurations were 
taken from the initial elastic equilibration run separated 
in time by at least 1000 time units. 



Ill INTERMEDIATE INCOHERENT 

SCATTERING FUNCTION AND SELF 

DIFFUSION CONSTANT 



v<(t) = Vi(t„) + / ds^(s) 



(4) 



The noise &(£) is Gaussian with zero mean and variance 



iW 



<rw(*0> = ^,i<w*(*-*') 



(5) 



,( a ) 



for the cartesian components Q , a = x,y,z. The 
stochastic process is implemented in the simulation by 



In this section we investigate time delayed correla- 
tions of a single tagged particle. In Fig. [T] we show for 
volume fractions 0.05 < r\ < 0.4 and for inelasticities 
e = 0.8,0.9,1.0 the incoherent intermediate scattering 
function 
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^incoh(q,i) = (^E eiq,(ri(t) " ri(0) 



(') 



— 


T| =0.05 E=O.S 





0=0.9 





£.= 1.0 





1\=0.1 


-- 


11=0.2 


— ■ 


r\=0.3 


-- 


\\=0.4 




300 



200 



100 



t ^ 




FIG. 1: Incoherent intermediate scattering function for sev- 
eral values of volume fraction r\ and coefficient of restitution e. 
All lines for N = 200000 where e = 0.8/0.9/1.0 are indicated 
with dashed/solid/dot-dashed lines respectively. r\ = 0.05 
corresponds to the left and r\ = 0.4 to the right lines. All er- 
ror bars are of the order of 10 -3 . Open diamonds and circles 
are for N = 10000, e = 0.8 and r\ = 0.05, 0.1 respectively. 



Since Fi nco h(q, t) is a measure of the correlation of par- 
ticle i at position ri(i) at time t and at position rj(0) at 
time t = 0, we find as expected that -F| nco h(q, t) decreases 
with increasing time. With decreasing densities the r, (i) 
and i"j(0) become more quickly uncorrelated and there- 
fore the decay is faster for smaller volume fractions. For 
the lowest densities the inelastic system can hardly be 
distinguished from the elastic case. For higher densities 
the relaxation is increasingly faster for the more inelas- 
tic systems. To quantify this effect we plot in Fig. [2] the 
relaxation time r when the incoherent intermediate scat- 
tering function has decayed to 1/e of its initial value, i.e. 
-F]ncoh(q, t) = 1/e. Clearly the elastic system shows the 
most rapid increase of relaxation time with density, even 
though the highest volume fraction (rj — 0.4) is still well 
below the critical value for the glass transition. The slow- 
ing down is weaker for the inelastic systems. However, 
the inelastic system also shows an increase by a factor of 
12 (e = 0.9) and 7 (e = 0.8). This indication for a pre- 
cursor of a glass transition even for the inelastic system 
is consistent with the higher density results of Kranz et 



27] 



al. [22j (theory and simulation) and of Reis et al. 
and Reyes et al. [1J] (experiment.) 

The intermediate incoherent scattering function 
-Pincoh(q, t), is often approximated by a Gaussian 



-Flncoh(q, t) 



-(Ar 2 (t)) 



assuming that the mean square displacement 
<Ar 2 (*)> = (iX>(t)-r i( 0)/ 



(8) 



(9) 



FIG. 2: Relaxation time of the incoherent scattering function 
as a function of volume fraction rj for several values of e for 
simulation runs with N = 200000 particles. 
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FIG. 3: Mean square displacement for e = 0.8 and various 
volume fractions. 



obeys Gaussian statistics. To test this hypothesis we first 
compute the mean square displacement (Ar 2 (i)). Fig. [3] 
shows the resulting (Ar 2 (t)) for e = 0.8 and several vol- 
ume fractions. One clearly observes a ballistic regime for 
small times with a crossover to diffusive behavior around 
t~ 1. 

The computed (Ar 2 (i)) is then substituted in the 
Gaussian approximation of Eq. ((5J) and compared to the 
full scattering function in Fig. [4] The Gaussian approxi- 
mation works very well for the densities under consider- 
ation, in particular for small q. 

We can also extract the self diffusion, 



D 



lim 

t— foo 



<Ar 2 (t)) 



6t 



(10) 




FIG. 4: Intermediate incoherent scattering function 
-Fincoh(q, t) using Eq. (solid lines) and for comparison the 
Gaussian approximation using Eqs. ((8]) and ((9| (dashed lines.) 




FIG. 6: Time-density superposition for the incoherent scat- 
tering function. 




FIG. 5: Self diffusion constant D as a function of volume frac- 
tion n. The filled symbols with lines for the eye are obtained 
via linear fits to the mean-square displacement for large times. 
Garzo results are corresponding to Eq. (2.10) of [23] and the 
Enskog result corresponds to Eq. (5) of [2ll |. 



via a linear fit to (Ar 2 (t)) at long times. The resulting 
D is plotted in Fig. [5] as a function of density (filled sym- 
bols) and compared with theoretical predictions (open 
symbols). As expected, the diffusion constant decreases 
strongly with density. Whereas the prediction of Enskog 
(see Eq. (5) of [2l|) is in excellent agreement for the elas- 
tic case (see inset), the prediction of Garzo [28J is very 
good for the inelastic case and r\ > 0.1. 

For the glass transition in elastic systems, one ob- 
serves dynamic scaling as the transition is approached. 
In other words, the scattering function does not depend 
separately on time and control parameter — either tem- 



perature or density — but only on the ratio t/r. We 
have tested this time-density superposition principle by 
plotting F inco h(q, t/r) for five volume fractions in Fig. [6] 
Even though the volume fractions under consideration 
are far away from the critical value, the data collapse for 
1] > 0.1. 



IV INTERMEDIATE COHERENT SCATTERING 

FUNCTION AND LONGITUDINAL CURRENT 

CORRELATION 

Information about the collective dynamics and in par- 
ticular collective density fluctuations is contained in the 
intermediate coherent scattering function, defined by 



1 



n 



F(q,*) = (^E e ^ Mt) ^' (0) 



(11) 



In the hydrodynamic regime, i.e. small wave numbers, we 
expect to see sound modes. This expectation is indeed 
born out by the data with an example shown in Fig. [7] for 
volume fraction 77 = 0.1, restitution coefficient e = 0.9, 
and for several q values. We observe oscillations which 
are overdamped for large q. 



A detailed analysis of the coherent correlation in terms 
of damped sound waves and temperature fluctuations will 
be given in Sec. VI, using a hydrodynamic model dis- 
cussed in the next section. Here we consider the longitu- 
dinal current correlation to obtain an approximation to 
the sound velocity which we analyze in dependence on 
volume fraction and inelasticity. The correlation of the 
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FIG. 7: Intermediate coherent scattering function F(q,t). 




FIG. 8: Longitudinal current correlation as a function of time 
for several densities; inset: variation with q. 



longitudinal current is denned as 



/ 1 N 1 



•faCO-r^O) 



iV " r 



r 



^(q,*). 



(12) 



In Fig. [3] we show data for two values of restitution and 
volume fractions 0.05 < r\ < 0.2. For all parameters we 
observe well-defined oscillations which are more strongly 
damped for the more inelastic system. 

In Fig. [9] we plot C\ (q, u) the corresponding Fourier 
transform of the current correlation. 

The position of the maximum of C/(q, w) can be used 
to estimate the speed of sound. The peak position w max 
as a function of wave number q is shown in Fig. 1101 As 
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FIG. 9: Longitudinal current correlation as a function of an- 
gular frequency. 
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FIG. 10: Position of the maximum, a; max of longitudinal cur- 
rent correlation C;(q,aj). Solid lines and filled symbols are 
for e = 0.9 and dashed lines indicate linear fits with slopes as 
listed in Table [TTT1 The inset shows for comparison u ma x for 
e = 0.9 (solid lines) and for e = 0.8 (dash-dotted lines). 



shown in the inset, the peak position does not depend on 
e. For small wave numbers a linear dispersion is observed 
(dashed lines in Fig. ITU]) , while deviations from linear 
behavior for larger wave numbers are more pronounced 
for the denser systems. 



V FLUCTUATING HYDRODYNAMICS 

In this section we compute S(q, u) from fluctuating hy- 
drodynamics. Our presentation follows closely the work 
of van Noije et al. [19( , except that we take care to con- 
serve momentum at each instant of time, thereby avoid- 



ing a divergence of the static structure factor. 

The hydrodynamic equations for the number density 
n and flow velocity u are the same as for an elastic fluid. 
However the equation for the temperature differs due to 
the energy dissipation in collisions and the energy input 
due to driving: 

d t T = D T AT - -^Vu - r + m& + 6. (13) 
an 

Here we present results in d = 3 dimensions. The en- 
ergy dissipation due to collisions, T, is estimated as 

i — 2 
r = 2Tv co \\ 2 J with the collision frequency v co \\. The 

input of kinetic energy due to driving is given by to£q, 

p denotes the pressure and Dt the thermal diffusivity. 

We have ignored nonlinear terms involving the flow field 

because we will consider only linear hydrodynamics. In 

the stationary state the energy dissipation in collisions 

and the energy input due to driving balance on average: 



and 



To = m£ . 



(14) 



We expand in fluctuations around the stationary state: 
n = n + Sn,T = T + ST and T = T + <5I\ The colli- 
sion frequency should be proportional to the density, the 
pair correlation function at contact, Xi an d the thermal 
velocity: v co \\ ex nxT 1 ' 2 , hence linearization around the 
stationary state T yields: T ~ r (l + ^ + i^<5n+f^). 
Following van Noije et al. [l9(, we consider a hydro- 
dynamic description of a granular fluid based on conser- 
vation of particle number and momentum and the re- 
laxation of temperature to its stationary value, To. The 
transverse momentum decouples so that we are left with 
three equations for the fluctuating density Sn, the longi- 
tudinal flow velocity u(q, t) — q-u/g, and the fluctuating 
temperature 5T: 



d t 6n{q,t) 
d t u(<i,t) 



-iqn u(<i,t) 

iq 



rano 

,2 



—5n{q,t) 
on 



pT M 



(15) 
(16) 



- vi<ru(q.,t)+€i(ti,t) 



2 P o 



d t ST(q,t) = -D T ^ST(ti,t)-iq-p-u{ti,t) 



dno 



(5n(q, t) f 1 d X 
Xdn 



"o 



<5n(q, t) 



(17) 

3 <JT(q, t) 

2 To 



0(q,i), 



where Dt = -3^- with the heat conductivity k, and where 
v\ is the longitudinal viscosity. Fluctuating hydrodynam- 
ics for an elastic fluid (e = 1) is based on internal noise, 
£( n and 8 ln , consistent with the fluctuations-dissipation 
theorem. Here we consider a randomly driven system: 
the particles are kicked randomly, giving rise to exter- 
nal noise in the equation for the velocity as well as the 
temperature. The external contributions are 



1 



rM) = -»i(*)<5(r. 



n) 



(18) 



9 cx (r,t) 



2m ■s—^ 

d^ ^ Vl 



ZidMr-Tiit)), (19) 



with variance 



c-2 

so 



(r(q,*)r(-q,t')> = v^5(t - t')(i - M) (20) 

no 



and 



: (q,^ cx (-q,0) = V^e 5(t-t')(l-5^ o ). (21) 
dno 



Here we have taken care of global momentum conser- 
vation, as realized by our driving mechanisms involving 
pairs of particles. These terms occur only for q = and 
ensure that the driving force vanishes at zero wave num- 
ber. Including both types of noise, £i — £} n 4- £, cx and 
0\ = 0[ n + 6>j cx , as suggested by Noije et al. [1J|, one 
obtains 

<£l(q,£)£l(-q,t')) =V(§- + ^^) S(t-t')(l-S^ ) 



no mno 



and 



(0(q,t)0(-c L ,t'))=W 



( mT $ 
\ dn 



2kT V 

~^PnT 



(22) 



<5(i-t')(l-<5 q ,o). 
(23) 



To complement the above equations, we need an ex- 
pression for the pressure p in terms of the density and 
temperature. Since the driven granular gas is far from 
equilibrium, we cannot expect that a thermodynamic 
description including an equation of state should hold 
in general. Nevertheless for small to moderate inelas- 
ticities an equation of state has been found empirically 
(Eq. (17.29) in [10]): p 2 nT (1 + 2n X (l + e)). We use 
the Carnahan-Starling approximation (for d — 3) 



X 



(1-77/2) 
(l-?7) 3 ' 



(24) 



This leaves us with two unknown parameters in the hy- 
drodynamic description, namely the longitudinal viscos- 
ity, v\, and the thermal diffusivity, Dt- 

The linearized equations can be solved for the fre- 
quency and wave number dependent correlation func- 
tions, S(q,u>) and Ci(q,oj) = ^ T S(q,uj). Of particular 
interest is the pole structure in the complex w-plane, de- 
scribing damped sound modes and the decay of temper- 
ature fluctuations. The latter can be either diffusive or 
with a finite relaxation rate, depending on wave num- 
ber q. For DTq 2 -C fjr-, the thermal diffusivity can be 
ignored and we have poles at 

, .3r 

zi 



UJ K 



±cq ± i'yq 



(25) 
(26) 



The latter correspond to sound modes with sound veloc- 
ity 



y th 



where v 2 = ± ( &■) 

tn m \ on J r 



2po L no 9x 

3mno \ x dn 



and damping 



(27) 



27 = v\ + 



4poT 
3romno ^3 



1 no_dx 
X dn 



Pi) 



dT n 



(28) 



In the opposite limit Drq 2 3> fjf 1 , we recover ordinary 
hydrodynamics of an elastic fluid. The sound speed is 
given by the adiabatic value 



M 



y th 



drnTQUy 



(29) 



and the temperature decay is diffusive ujt — ±iDTq 2 -%- ■ 
In general, we expect to see a crossover, when Dxq 2 = 
|j£. In order to estimate q c , we use the Enskog values 
for the collision frequency in three dimensions and the 
thermal diffusivity: 



Veda = we = 47rxn a' 



Dn 



75 



d{\ + e)(49 — 33e)n a 2 x V nm 



To 



(30) 
(31) 



These yield for the crossover wave number 



2 2 _ 6(l-^)(l + £ )(49-336)xW ,„„, 



Numerical values of estimated q c a for the simulated vol- 
ume fractions r\ and inelasticities e are given in Table U 
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e 


q c a 


0.05 


0.8 


0.21 


0.05 


0.9 


0.14 


0.1 


0.8 


0.48 


0.1 


0.9 


0.33 


0.2 


0.8 


1.29 


0.2 


0.9 


0.89 



TABLE I: Estimates for q c a using Eq. ((32]l . 



For the case of Drq ~ |4£ it is necessary to use the 
more general solution for the dynamic structure factor 



S(q,u)) = n q 2 



lo 2 + (370WE + D T q 2 Y 



§1 , qujjhr 
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where we have used i/ co \\ = we, the abbreviation 70 = 2 J and where 
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<4^t<T 



(34) 



VI COHERENT SCATTERING FUNCTION AND 
TRANSPORT COEFFICIENTS 



According to the g c estimates given in Table HI we see 
that our data arc neither clearly in the hydrodynamic 
regime nor in the inelastic regime, but in general the two 
relaxation terms in the equation for the temperature are 
comparable in magnitude. Hence, we fit the simulation 
results of the dynamic structure factor to the full ex- 



pression for S(q,u) as given in Eqs. (|3"3")l and (f3~4"|) . We 
allow for two fit parameters, Dt and vi, with all other 
parameters determined by the approximate equation of 
state. The best fits (solid line) are shown in Figs. QT] 
- 1141 in comparison with the simulation data (symbols) 
for S(q,cu). We find excellent agreement not only for 
very small <?, for which we would expect best agreement 
with the hydrodynamic equations, but also for q < 1.0. 
Both features, the shoulder due to the sound wave as 
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FIG. 11: Dynamic structure factor S(q,uj) for rj = 0.05, e = 
0.8 and q — 0.2 - 0.5. Symbols indicate simulation results 
obtained via Fourier transform of F(q, t) and lines indicate 
fits with Eq. (J33J). 



FIG. 13: S(q,co) for rj = 0.2, e = 0.9 and q = 0.5 - 0.7. 




FIG. 12: S(q,uj) for rj = 0.05, e = 0.8 and q = 1.0 - 3.0. 



well as the damping, are quantitatively in agreement with 
Eqs. (|33|) and |34|) . Similarly, we find very good agree- 
ment for the r\ = 0.1 results. 

The corresponding best fit parameters are the trans- 
port coefficients Dt and t'l which are shown graphically 
in Figs. [TBI fe 1161 The fits require g-dependent transport 
coefficients because we consider wave numbers outside 
the hydrodynamic regime. It is difficult to estimate the 
hydrodynamic regime, but we need at least q < q c (see 
Table El, corresponding to D T q 2 < f^. For rj = 0.2, 
we are able to reach this regime and indeed find that Dt 
and v\ are approximately independent of q. For rj = 0.05 
even the smallest q values are not in the hydrodynamic 
regime yet, and for rj = 0.1 the smallest wave numbers 
are in the crossover regime. As far as temperature fluctu- 




FIG. 14: S(q,uj) for r/ = 0.2, e = 0.9 and q = 0.8 - 1.5. 



ations are concerned, the diffusive regime is restricted to 
larger wave numbers Dxq 2 > |5r-, so that D T can only 
be extracted from an intermediate range of q- values, such 
that q > q c but q still small enough to ignore higher or- 
der terms in q. Again, for rj = 0.2 this seems possible, 
whereas for rj = 0.05 our data are not sufficient. 

Tabulated in Table [IT] is a quantitative comparison of 
the fit results for small q with the theoretical predictions 
for D T = %; and v x = \ ( j&tes^ + () , where r; shoar 

and C, are shear and bulk viscosity respectively. For the 
comparison with Brilliantov et al. we use Eqs. (20.13) 
& (20.30) of Ref. llOl f or ?7 s hcar and k respectively and C, 
of Eq. (32) of RefTlli. For the predictions of Dufty ct 
al. we used Eqs. (29), (30) and (32) of Ref. [H] and for 
the predictions of Garzo et al. we used Eqs. (Bl), (2.2), 
(3.8), (2.3) and (3.9) of Ref. 0. 

We find best agreement between the simulation results 
for the smallest q and the predictions of Dufty et al. [29( 
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FIG. 15: Thermal diffusivity Dr(q) obtained via fits to 
S(q,u>) with Eq. (|33[1 . The arrows indicate the q c values from 
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FIG. 16: Longitudinal viscosity v\(q) obtained via fits to 
S(q,u>) with Eq. (|33[1 . The arrows indicate the q c values from 
Tabled 



77 = 0.05 




£ = 0.8 


£ = 0.9 




D T 


V\ 


D T 


V\ 


Fit Results: q = 0.2 


4.72 


2.55 


4.63 


3.23 


g = 0.3 


3.34 


1.69 


3.45 


1.81 


g = 0.4 


2.69 


1.39 


3.00 


1.52 


Brilliantov et al. [10] 


3.19 


2.26 


3.54 


2.25 


Dufty et al. [29] 


4.71 


2.82 


4.07 


2.77 


Garzo et al. [28] 


5.62 


2.78 


5.06 


2.67 


77 = 0.1 




£ = 0.8 


£ = 0.9 




D T 


V\ 


D T 


V\ 


Fit Results: 5 = 0.3 


2.23 


1.20 


2.67 


1.70 


g = 0.4 


2.25 


1.02 


2.42 


1.25 


5 = 0.5 


2.15 


1.07 


2.33 


1.20 


Brilliantov et al. [10] 


1.39 


1.13 


1.55 


1.13 


Dufty et al. [29] 


2.67 


1.69 


2.42 


1.71 


Garzo et al. [28] 


2.81 


1.53 


2.53 


1.48 



?7 = 0.2 




£ = 0.8 


£ = 0.9 




D T 


V\ 


D T 


V\ 


Fit Results: 5 = 0.5 


1.95 


1.02 


2.22 


1.10 


5 = 0.6 


2.09 


1.02 


2.32 


1.24 


5 = 0.7 


2.20 


1.23 


2.27 


1.30 


Brilliantov et al. [10] 


0.52 


0.83 


0.57 


0.85 


Dufty et al. [29] 


2.03 


1.63 


2.01 


1.72 


Garzo et al. [28] 


1.40 


1.15 


1.26 


1.15 



TABLE II: Comparison of theoretical predictions and fit re- 
sults for Dt and V\. 





via Wmax(g) 


DW « i£ 


DW » §£ 


c 


= 0.8 


77 = 0.05 


1.58 


0.73 


1.55 


77 = 0.1 


1.81 


0.90 


1.87 


77 = 0.2 


2.57 


1.37 


2.79 


c 


= 0.9 


77 = 0.05 


1.62 


0.74 


1.56 


77 = 0.1 


1.89 


0.92 


1.90 


77 = 0.2 


2.66 


1.40 


2.86 



TABLE III: The speed of sound, c, 
the simulation results for oj max (5) 1 



determined via the slope of 
see Fig. llO|l compared with 



the predicted values of Eq. 1)270 in the case of Drq <C 
and with Eq. ([29j) in the case of D T q 2 > fp-. 



:<r„ 

2T„ 



and fairly good agreement with the predictions of Garzo 
et al. [28j . 

Finally, we compare the speed of sound as obtained 
from the maximum of the current correlation with the 
predictions from the hydrodynamic theory in either the 
inelastic regime (see Eq. ([27]) ) or the diffusive regime 
(see Eq. p9| L We find very good agreement (see Ta- 
ble [nlj of the simulation results with Eq. (|29|) . implying 
Drq 2 ~> ffJ 1 and adiabatic sound propagation. However 
one should keep in mind that our procedure to extract 
the sound velocity from the maximum of the current cor- 
relation yields only an estimate of the sound velocity. 



CONCLUSIONS AND OUTLOOK 

We have investigated a homogeneously driven granu- 
lar fluid of hard spheres at intermediate volume fractions 
0.05 < 77 < 0.4 and for constant normal restitution co- 
efficients 0.8 < e < 1.0. Using event-driven simulations 
we have determined time-delayed correlation functions in 
the stationary state. 

We find for the incoherent intermediate scatter- 
ing function that it follows time-density superposition 
and that it is well approximated by the Gaussian 

i ? incoh(q, t) = e~z~\ Ar2(t ^ , where (Ar 2 (i)) is the mean 



10 



square displacement. The decay time of Fi nco h(q, t) is 
rapidly increasing with increasing 77, giving rise to a cor- 
responding decrease of the diffusion constant. This pre- 
cursor of a glass transition, which occurs at significantly 
larger 77, is thus present not only in the elastic fluid but 
also in the inelastic case consistent with previous results 
at larger densities [13, [H, Hi | • 

We also determine the coherent intermediate scatter- 
ing function F(q,t), the longitudinal current correlation 
function C\{q,t) 1 and their Fourier transforms S(q,ui), 
C\(q,ui). Because we are interested in the long term dy- 
namics we have simulated comparatively small systems 
of N — 10000 particles and averaged over 100 indepen- 
dent simulation runs. We observe sound waves in the 
form of oscillations in F(q, t) and estimate the sound ve- 
locity from the peak of C\(q,t). For a quantitative com- 
parison with the predictions of generalized fluctuating 
hydrodynamics, we use the linear hydrodynamic equa- 
tions of Noije et al. [19J] and fit the solutions thereof to 
the simulation results for S(q,ui). Depending on wave 
number and inelasticity the temperature fluctuations are 
predicted to be governed by inelastic collisions or diffu- 
sion [19|, |30( • Our results are consistent with being in the 
"standard regime" [l9( in which the speed of sound is the 
same as for elastic particles (see Table [TTTj) and the damp- 
ing of the sound wave depends on inelasticity. The most 
accurate fits were obtained assuming generalized hydro- 
dynamic equations which account for both temperature 
diffusion as well as dissipation due to inelastic collisions 
(Eq. (|33| ). The resulting transport coefficients D, Dt 
and v\ compare well with the predictions of Dufty et al. 
[29J | and Garzo et al. [H]. 



We conclude that the time delayed correlations of a 
fluid of inelastically colliding particles are well described 
by generalized hydrodynamics. It would be interesting to 
extend our study in several directions. First, one would 
like to see still smaller q, requiring significantly larger 
systems (and yet also many independent runs for suffi- 
cient statistics). Second, it would be interesting to go to 
higher density and study sound propagation as the glass 
transition is approached. Finally, time- or frequency- 
dependent response functions are largely unexplored. 

K.V.L. thanks the Institute of Theoretical Physics, 
University of Gottingen, for financial support and hos- 
pitality. We thank Till Kranz for many interesting dis- 
cussions. 
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